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Abstract 

The dyadic Green's function for a current embedded in 
a grounded dielectric slab is used to analyze microstrip 
lines at millimeter wave frequencies. The dyadic Green's 
function accounts accurately for fringing fields and 
dielectric cover over the microstrip line. Using 
Rumsey's reaction concept, an expression for the 
characteristic impedance is obtained. The numerical 
results are compared with results, reported by others. 
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Introduction 


In recent years there has been a great deal of interest 
in analyzing microstrip lines at millimeter wave 
frequencies due to use of these lines in monolithic phased 
array systems. 

Microstrip lines were initially analyzed using a 
quasi-static approach and later by a waveguide model to 
study dispersion characteristics [4]. However, these 
models do not take into account losses due to radiation and 
surface wave excitation. These losses can be significant 
at high frequencies and some attempts [1-3] have been made 
to._aec.Qunt for these losses. In reference 1, the current 
distribution on the microstrip discontinuity is determined 
by solving the electric field integral equation. From a 
knowledge of the current distribution, characteristics such 
as impedance and guide wavelength can then be determined. 
The microstrip line embedded in a multilayer dielectric has 
also been analyzed using a generalized spectral domain 
Green's function [3]. 

In this paper, the dyadic Green's function for a 
current source embedded in a grounded dielectric slab is 
used to determine the field due to the microstrip line 
current. Using Rumsey's reaction concept [5], the 
characteristic impedance of the microstrip line embedded in 
a dielectric slab is then determined. The numerical 
results obtained using the present method are compared with 
earlier published data. 
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Symbols 
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<x 
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Ay 
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Total thickness of dielectric. 

Electric field vector. 

Incident electric field. 

Frequency. 

Transverse distribution function. 

Complex amplitude of n th pulse. 

Total input current. 

Surface electric current density vector. 

. 

x-component of J(x,y). 

Bidimensional Fourier transform of J^(x,y). 

Surface current distribution on microstrip. 

Surface current distribution on microstrip 
induced by x-polarized plane wave. 

Wave propagation constant in dielectric, k Qv / c 

Wave propagation constant in free space. 

Dominant mode propagation constant. 

Fourier transform variable with respect to x. 

Fourier transform variable with respect to y. 

Number of pulses in y-direction. 

Width of microstrip line. 

Unit vector along x-,y-,z-axis respectively. 
-W/2 + ( 2n-l ) Ay/2 . 

Characteristic impedance of microstrip line. 
Microstrip location from ground plane. 
Imaginary part of K . 

e 

Real part of K . 

e 

W/N. 

Relative permittivity of dielectric substrate. 
Wavelength in free space. 

Characteristic impedance of free space. 
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Theory 


The geometry of a microstrip line embedded in a 
grounded dielectric slab is shown in figure 1 along with 
the notation to be used. In order to determine the 
characteristic impedance, consider a microstrip line being 
excited by a current source backed by a perfect electric 
conductor at the x=0 plane, as shown in the current source 
looks into a semi-infinite transmission line, the 
characteristic impedance of the line would be the same as 
the input impedance seen by the source. Assuming the strip 
width, W, is much smaller than the operating wavelength, 
the y-component of the strip current can be neglected. The 
surface current, J(x,y), on the z=z' plane may then be 
represented by 

J(x,y) = xf(y) exp(-jKjx|) -oo^x^co (1) 

where K =/3+ ja is the dominant mode propagation constant 

e 

along the strip. The range of x is extended to due to 

the image of the strip behind the conducting plane at x=0 
and the magnitude of x is used to indicate propagation away 
from the source at x=0. The transverse distribution, f(y), 
in equation 1 may be assumed to be of the form 

N 

f(y) * £ I P n (y) (2) 

n=l 

where 

p„(y) = !/ A y for (y n - -j^) s y - (y n + -j^) 

= 0 otherwise. 

From the continuity condition, the distribution f(y) 
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roust satisfy 


W/2 

| f (y) dy = i 0 (3) 

-W/ 2 

The pulse amplitudes, I , in equation 2 are determined, 

n 

as explained later, by using the integral equation method 
in conjunction with the method of moments. 

If E(x,y, z) is the electric field due to the strip 
current, J(x,y), then the input impedance seen by the 
source is obtained from 


Z 

In 


Z 


0 



W/2 oo 

J | J(x , y) • E(x,y,z') dx dy 

-W/2 0 


( 4 ) 


This approach thus differs from the solution to the 
multilayer transmission line problem [3], where the 
characteristic impedance is obtained after dividing the 
average voltage by the total current on the strip. 
Furthermore, by taking into account the actual current 
distribution in the y-direction, it is expected that the 
present approach would give more accurate results. Using 
the dyadic Green's function for a current embedded in a 
grounded dielectric slab [6], the x-component of the 
electric field E(x,y,z) is obtained as 


x • E(x, y , z' ) 


— j rj k z' 
J 0 0 

( 2 7T ) 2 k 2 


00 00 

J j Q ( k „' k y) •M’VS) 

— 00 — 00 


where, 


• exp( jk x x+jk y y) dk^ dk y (5) 
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Q(V k y ) = (* a - k X(V k y ) + J k X(V k y ) 


( 6 ) 


M k * k ) = 

1 v x y ' 


sin (k r z ' ) 


(M') 


f kjCos (kj (d-z' )) + jk ij sin(k i (d-z'))’| 


k i cos(k i d) + jk Ji sin(k ] d) 
-k 2 s in (k f z ' ) 


S' (k ,k ) = 

2 x y 1 c ^ kj cos (k i d) + jk j sin(k i d) 


sin (kj z ' ) ’ 

f (Gr_1) ] 

l (V) J 

1 k J cos(k i d) + jk ii sin(k i d) J 


(7a) 


(7b) 


y k 2 -k 2 -k 2 


X y 


l -j/k 2 +k 2 -k 2 


x y 


(k 2 +k 2 ) * k 2 

v x y ' 


(k 2 +k 2 ) > k 2 

v x y J 


II 


i 2 i 2 i 2 

k -k -k 

0 x y 


-j / k ? +k 2 -k 2 

J V x y 0 


(k^+k^) ^ k^ 

' x y' 0 

(k 2 +k 2 ) > k 2 

v x yf 0 


In equation 5, tj is the free space characteristic 
impedance/ 1?^ and are derived in reference 6/ and 
J^(k^,k y ) is given by 

W/2 co 

J x( k x ,k y) = J J f ( y > ex P("i K e l x l“3 k x x ‘j k y y) dx d Y (8) 

W/2 -oo 

From equations 2 and 8, and assuming a to be a finite 
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negative value, the integral in equation 8 may be evaluated 
in closed form as 



N 


I 

n=l 



sin (k y Ay/2) 

( k /y/ 2 ) 



1 

k +K 

x e 


exp(-jk y yj 


( 9 ) 


Substituting equation 5 into equation 4 and after some 
mathematical manipulations, the characteristic impedance is 
obtained as 


Z = 
o 


, CO CO 

j 7) k z N 1 1 r r 

£ Hr I J 

n=l I J J 

0 —oo -co 


( 2tt ) 2 k 2 


sin (k y Ay/2) 


L ( k v A y/ 2 ) 


k -K £ +K 

x e x e 


Q ( k x' k ,) 


J f-k ,-k ) exp(-jk y ) 

^ X V x y ; v y n ' 


[ dk dk 


x y 


( 10 ) 


where 


W / 2 oo 


J x (-k x ,-k y ) = J | f (y) exp(-jK p x + jk x x + jk y y) dx dy 


-w/ 2 0 


N 

E J 1 . 

m= 1 


sin (k y Ay/2) 


(k y iy/2) 


eK P( 


k -K 

X e 


( 11 ) 


Substituting equation 11 into equation 10 gives 
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z 


o 


-jTj k z' N Nil 

0 0 ^ „ n m 

? ? L. L, P 

(2rc) k n= 1 m= 1 I 


03 CO 

r r i 


sin (k y Ay/2) " 

2 

1 1 

j j 
00 “00 

l 

. ( k /y/ 2 ) . 


(k -K ) 2 k 2 -K 2 

L v x e J x 


• Q ( k x ' k y ) exp(j k y ( yni -y n )) 


dk dk 

x y 


( 12 ) 


Making a change of variables such that: 


C - 




K 



o 


equation 12 is rewritten as 


Z 

o 


-jT) k z' N N 

— i i 

(2tt) k n=l ra= 1 


I I 

n m 



00 00 , 

1 1 

00 “CO k 

"sin (£ k o Ay/2)' 

2 

r i 

1 

. (s k 0 A y/ 2 ) . 


i 

OJ 

w2 2 

C -v J 


• Q( C,€) e x p(je k 0 (y m -y n )) 


dC d£ 


(13) 


where, Q(C,£) ■ Q(k Q C,k o C) . 

Note that the integrand in equation 13 has poles at 
of order two and simple poles at C=±i/. The integration with 
respect to £ can be evaluated using the residue theorem. In 
order to evaluate the integration with respect to £, 
consider the following contour integration, 
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(14) 


J{ d} 

where C is the real part of % , and C is the closed contour 
to be selected. The function <?(?,£) has branch cuts at 

* - * J i - ? 2 

and poles corresponding to the roots of the equations: 

k cosfk d) + jk sin(k d) = 0 
c k cos(k d) + jk sin(k d) = 0 

These poles correspond to the surface wave inodes that would 
be launched by the serai-infinite microstrip line. These 
surface wave poles and the branch cuts may be plotted in 
the complex of-plane as shown in figure 3. The contour C 
should be selected such that the propagation of the surface 
waves in the x>0 direction is assured. For this condition, 
the surface wave poles on the positive real axis must be 
included in the contour C. The required contour C is shown 
in figure 3. The integration in equation 14 can be written 
as 

J = J+ J+ J= -2nj { R + R 2 + R 3 ) (15) 

C C C C 

t 2 3 

where, 

R = E Residues of Q at surface wave poles 
R 2 = Residue of Q/(q-v) at %-v 
R 3 = Residue of Q/(cj 2 -i 2 2 ) at %=v . 

In equation 15, the integration along C 2 is zero because of 
the radiation condition, the integration along the branch 
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cut, C 3 , represents a radiation field, and the integration 
along C is the required integration to be evaluated with 
respect to C, . Since the radiation field due to the 
microstrip line would be very small, the integration along 
C 3 will be neglected. For the present case, the 
contribution due to surface wave inodes will also be 
neglected. The integration with respect to C becomes 



(^ 7 ) } 0<c ' e > dc 


= — 2tt j 






<?(C,Q 

2C 


t=v 


} 


(16) 


Substituting equation 16 into equation 13, the expression 

for Z reduces to 
0 


00 


7| k z' N N II 

0 0 j-. n m 

* 


'sin (k o ^Ay/2)' 

2nk 2 n=l m=l I 2 

0 
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. (k 0 CAy/2) _ 


0 


23Q(C,C) 


Q(CfC)’ 




cos (e(m-n)k Ay) i d$ 


(17) 


From the characteristic equation (27), it can be seen that 
the term 


N N 


E E 

n=l m=l 



0K,£) 

C 




'sin (k Q CAy/2) 


(k 0 ?Ay/2) 


cos (£(m-n)k o Ay) 


d£ 


in equation 17 is always^ zero. 
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where J Q (*) denotes the Bessel function of the first kind 
of zero order. 
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Evaluation of K and I 

e n 

In order to determine K , consider the infinite 

e 

microstrip line as shown in figure 1. The surface current 
distribution on the strip may be assumed to be of the form 

J t (x,y) = * J xt ( x ' y ) = x exp(-jK e x) — X— oo (21) 

where f(y) is as given in equation 2. 

Using equation 8, the Fourier transform of J t (x,y) is 
obtained as 



where 5(*) is a delta function and is the result of 
assuming a traveling wave current in the x-direction on the 
strip as denoted by the exponential factor in equation 21. 
The x-component of the electric field radiated by J 

Xt 

is obtained by replacing J (k , k ) with J (k ,k ) in 

xx y xt x y 

equation 5. Then by equating the x-component of the 
electric field to zero on the surface of the strip, since 
the total field must be zero and there is no incident 
field, an electric field integral equation is obtained 
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-It} k z' N 

0 = — — I I 2n 

( 2tt ) 2 k 2 n=l n 


CO CO , 

r r J 

sin (k Ay/2) 

J J 

CO —00 k 

( k y Ay/2) J 


Qfk ,k ) s(k +K ) 

v x y 7 v x e 7 


exp(-jk v n ) exp(jkx+jk y) 


dk dk 


( 23 ) 


By setting x=0 in equation 23 and selecting P m (y) as a 
testing function and multiplying both sides of equation 23 
by P (y) yields 


-Itj k z' N 

0 = — - — E I 2tt 

( 2tt ) 2 k 2 n=l n 


00 00 , 

r r 

sin (k y Ay/2) 

J J 1 
00—00 k 

_ ( k y Ay/2) 


Q( k ,' k y ) «(V K .) 


P.(y) exp(-jk y yj exp(jk y) 


dk dk 

x y 


(24) 


Integrating both sides of equation 24 with respect to y and 
rearranging gives 


0 = 


-Itj k z' 
_ 0 0 

( 2tt ) 2 k 2 


N 

I 
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V" I I 


sin (k y Ay/2) 
(k y Ay/2) 


• Q( k x ' k y ) 6 ( k x +K J ex pH k y yJ 


y+Ay/2 

m r 

p (y) ex p(i k y) d y 

- y m -Ay/2 


dk dk 

x y 


(25) 
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Using the property, 

00 

| f(x) <5(x+x' ) dx = f(-x') 

“ 00 

of a delta function, equation 25 can be written in the 
form 


0 = 


— j t? k z' N 
J 0 0 £ 

2rr k 2 n=l 


I 

n 



k =-K 


:in (k y Ay/2) 
(k y Ay/2) 


i 2 


exp(- jk (y n -yj) 


dk 


(26) 


Since Q(k x ,k y ) =Q(k x ,-k y ) , equation 26 may be rewritten as 

m 

■jT] z' N 


0 = 


iff Cj r 

ITT- 

r n=l 




sin (Ck 0 Ay/2) 
. (Ck Ay/2) 


cos (Ck Q Ay(m-n)) 


dC (27) 


where £=k /k^, ^=k x /k Q , v-K^/k^, and Q(C f C)=Q (^ 0 C / • 

By selecting m=l,2,3 N, equation 27 may be written 
in matrix form as 
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2 t 
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22 


z 


. 32 


z 


N2 


z 

13 

z 

IN 


1 

M 

z 

23 

z 

2N 


I 

2 

z • • • 
. 33 

z 

. 3N 


I 

3 

z 

N3 

z 

NN J 


I 

L N 
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where the components of the square matrix are given by 


Z 

mn 


" r f , , 


sin (£k Q Ay/2)' 

J Q«,C) 
0 ' 


. (Sl< 0 Ay/2) . 


* cos (£k Q Ay(m-n) ) > 

; 


dC 


(29) 


For a nontrivial solution of equation 28, the determinant 
of the square matrix in equation 28 must be set to zero, 
therefore, 
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1 1 
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z 
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z 

32 
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. 3N 

7 

N 1 

z 

N2 

7 

N3 

7 

NN 


Equation 30 is solved for K e in a numerical procedure 
by selecting N and then starting with an estimate of K 
and iterating K until equation 30 is satisfied to the 

e 

desired accuracy. 

In order to determine the current distribution in the 
transverse direction, it is assumed that the infinite 
microstrip line shown in figure 1 is excited by an 
x-polarized plane wave. The surface current induced on the 
strip may be assumed to be of the form 


J(x,y) = x J xu (x,y) 


N 

l I P (y) x 

- n n 


(31) 


If E (x,y,z) is the x-component of the scattered 

XU J 

electric field due to J^(x,y), and E x (x,y,z) is the 
incident electric field, the boundary condition of zero 
tangential electric field on the conducting strip gives 
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( 32 ) 


-e'(x, y , z' ) = E (x,y, z' ) 

X 

E (x,y , z' ) can be obtained from equations 5,6 and 8, 

xu' 

Equation 32 then becomes 


-E*(x,y, z' ) 


- j T) k z ' 
J 0 0 

(2n) 2 k 2 
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l 

n=l 


I 
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co 


I 



k =0 
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sin (k Ay/2) 


(k Ay/2) 


exp(-jk yj 


dk 


(33) 


Equation 33 is an integral equation with 1^ as the 
unknown complex values to be determined. Equation 33 can 
be solved by selecting P (y) as a testing function and 

J m 

applying the method of moments to yield 
N 

V I z = - Ay e‘(z') for m=l , 2 , 3 , ■ • ■ ,N (34) 

^ - n mn 0 


where E^(z') is the x-polarized plane-wave electric field 
incident on the strip from the z-direction (normal 
incidence) and Z is given by 

— • 7 mn 


z 

mn 


- jf) k 2 z' Ay 
J 0 0 J 

jjY — 

r 


oo 
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<.’(05) 


C=o 


sin (5k o iy/2) 

(?Vy/2) 


cos (^(m-n)k Q Ay) | d£ 


( 35 ) 


Equation 34 can be written in matrix notation as 
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( 36 ) 


and the complex values of 1^ are obtained by matrix 
inversion . 
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Numerical Results 


In order to obtain numerical results, the semi-infinite 
integrals appearing in equations 17, 18, 20, 29 and 35 are 
evaluated using Gauss Quadrature integration techniques and 
the upper limit on the integration is truncated to a large 
finite value. From numerical results, it was found that 
the upper limit on the integration depended upon frequency 
and could vary from t~200 (for high frequencies) to t~1000 
(for frequencies below 5 GHz), Because of the complexity 
of the expression for Q(k ,k ), the derivative, 9Q/3k or 
dQ/d is evaluated numerically. 

In order to obtain meaningful numerical results, it is 
desirable to apply convergence tests to 1^, and Z Q . The 
matrix equation is solved for I ( for various values of N 
and plotted in figure 4. It is observed that a stable 
solution for I is obtained when N^25. Figure 5 shows the 

n 

variation of K and Z as a function of N. It is noted 

e 0 

that N— 2 5 also gives converged values for these quantities. 

To verify the present formulation, the characteristic 
impedance and the guide wavelength of a microstrip line 
with W/z' *1.0, d=z'=-.04 cm, and e =10, is computed using 
equations 13 and 22 for N=25. These results are presented 
in figure 6 as a function of frequency along with results 
reported earlier [5J. Using the present method, Z q and 
are also computed for uniform current distribution and 
nonuniform (see eq.14) current distributions and presented 
in figure 6. The results obtained by the present method 
with uniform distribution are in excellent agreement with 
the results reported earlier [4]. However, Z q and K g 
calculated by taking into account the actual current 
distribution (as obtained from equation 36) is believed to 
give more accurate results than the results obtained with 
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uniform current distribution. It is also observed from 
curves C & D of figure 6, that the assumed nonuniform 
distribution is close to the actual calculated current 
distribution and may be used instead of the pulse 
distribution in order to save computation time. To compare 
the present method with the approach in reference 3, the 
characteristic impedance of a microstrip line with d=z and 
d=4z' e =2.53, f =3GHz is computed as a function of strip 

r o 

width, W/z' , and presented in figures 7 and 8. It is noted 
that there is a small discrepancy between the two results 
in figure 8. This may be attributed to the different type 
of current distribution assumed in the two cases. With an 
assumption of uniform current distribution in the 
y-direction, the results obtained by the present method are 
also presented in figure 7. There is good agreement 
between the results obtained by the two methods. Figure 8 
shows the characteristic impedance of a microstrip line 
covered with a dielectric slab. 

The results obtained by the present method are also 
compared in figure 9 with data taken from reference 1. For 
frequencies greater than 42 GHz, the dielectric thickness 
becomes greater than a quarter-wavelength and the current 
distribution in the y-direction may not be of the assumed 
form. Furthermore, for thicknesses greater than a 
quarter-wavelength, the microstrip line would act as a 
radiator rather than a transmission line. 

It is interesting to study the effect of distribution 
in the y-direction on the guide wavelength and the 
characteristic impedance of a microstrip line. It is 
expected that for a narrow strip, the distribution may be 
closer to the assumed nonuniform function (equation 19); 
whereas, for a wider strip, the distribution is expected to 
be closer to uniform, except near the strip edges. To 
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verify this expectation, the guide wavelength and the 
characteristic impedance of a microstrip line with 
d=z'=-.05 wavelengths, e=10, and the strip width in the 
range of 0.025 to 0.330 wavelengths were computed and 
presented in figure 10. The guide wavelength computed with 
the assumed nonuniform distribution and pulse distribution 
are almost the same for a narrow strip (W/A o <0 . 05) . For 
wider strips (W/A q > 0.27) the guide wavelength computed 
using the uniform distribution and the pulse distribution 
approach the same values. From figure 10, it is observed 
that the characteristic impedance does not significantly 
depend upon the transverse distribution of the current. 

For microwave integrated circuit design, alumina 
substrate with dielectric constants in the range of 10 to 
12.8 are often used. In order to design a microstrip line 
on these substrates, an accurate knowledge of the guide 
wavelength and the characteristic impedance is essential. 

In figures 11 and 12, the guide wavelength and the 
characteristic impedance are presented as a function of 
frequency for various design parameters. Calculations for 
other parameters may be obtained with the aid of the 
computer code given in the appendix. 
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Conclusion 


Using the dyadic Green's function and Rumsey's reaction 
concept, an expression for the characteristic impedance of 
a microstrip line embedded in a grounded dielectric slab 
has been derived. The characteristic impedance determined 
by the present method is found to be in excellent agreement 
with earlier published results. Numerical results for a 
microstrip line on alumina substrate is presented at 
millimeter wave frequencies. A computer code listing is 
included for calculating the characteristic impedance and 
propagation constant for microstrip lines. 
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Figure 2. Current source exciting microstrip line. 
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Figure 3. Complex 7 -plane with surface wave poles and branch cuts. 
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Number of Pulses 


Figure 5 


Variation of characteristic Impedance and normalized guide 
wavelength of mlcrostrlp line as a function of number of 
transverse pulses, (c = 10, W^z', z'*d™0.04 cm, f*60 GHz). 





Wavelength / Guide wavelength 



Frequency (GHz) 


Figure 6. 


Characteristic Impedance 
of mlcrostrlp line as a 


( 


uniform, - - - 


and normalized guide wavelength 
function of frequency, 
nonuniform, N=25 pulses). 


Figure 7. Characteristic 1 Impedance of covered mlcrostrlp line, 
(d * 0.64 cm, z'- 0. 16 cm, c - 2. 53, f ■ 3 GHzl 






Characteristic Impedance (o 



W/z' 

Figure 8. Characteristic impedance of microstrlp line. 
(z'» d - 0.16 cm, e r « 2.53, f = 3 GHz). 
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Characteristic Impedance (ohms) 



Wavelength / Guide wavelength 





Figure 10. 


Characteristic Impedance and normalized guide wavelength 
of mlcrostr ip line as a function of strip width. 

(e ■ 10, z'* d • 0.05 

(•^ uniform, •- nonuniform, — N=25 pulses). 
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Characteristic Impedance (ohms) 



Wavelength / Guide wavelength 





Frequency (GHz) 


Figure lib. Characteristic Impedance of mlcrostrlp line as a 

function of frequency versus strlp-wldth/strlp-helght 
(e * 9.6, z'* d * 0.06 cm). 





Wavelength / Guide wavelength 





Figure 12a. 


Normalized guide wavelength of mlcrostrip line as a 
function of frequency versus strlp-wldth/strip height, 
(c « 12.8, z'* d » 0.06 cm). 

r 
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Frequency (GHz) 


Figure 12b. Characteristic Impedance of mlcrostrlp line as a 

function of frequency versus strlp-wldth/strlp-helght 
(c » 12.8, z'» d - 0.06 cm). 
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Appendix 

Instruction for using computer code 


The main program MSTRIP.FOR calculates the characteristic 
Impedance, Z^, and the dominant mode propagation constant, K^, of 
a mlcrostrlp line embedded In a grounded dielectric slab. 


The main program must be linked to the following subroutines: 

ZTEMP. FOR 
SRMN . FOR 
HATINV. FOR 
MATROW . FOR 
RSTRIP1. FOR 
RSTRIP. FOR 
MZO.FOR 
DFKX. FOR 

INPUT DATA : 

NU = Integer selecting the type of current distribution. 

1 for pulse distribution. 

2 for uniform distribution. 

3 for nonuniform distribution. 


MM 

ED 

EZP 

E 

El 

WSTRIP 

WMAX 


= Number of pulses assumed In the pulse distribution, 
(maximum number = 242). 

a. Dielectric slab thickness In free space wavelengths. 

* Height of strip above ground plane In free space 
wavelengths. 

= Dielectric constant of slab. 

= Loss tangent of dielectric slab material. 

= Minimum mlcrostrlp width in free space wavelengths. 

= Maximum mlcrostrlp width In free space wavelengths. 


OUTPUT DATA : 

For pulse distribution, the transverse current distribution Is 
output to file F0R010.DAT. 

The characteristic impedance and dominant mode propagation 
constant as a function of mlcrostrlp width are output to file 
F0R013.DAT. 
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The current distribution is assumed in the form; 


J x (x,y) - f(y) exp(-JK^x) 

J y (x,y) = 0. 

For the above current distribution, the program calculates the 
dominant mode propagation constant r and the characteristic 
impedance for one of three cases: 


f(y) - io/w 


(uniform) 



(nonunlf orm) 


( pulses ) . 


For the distributions, the program calculates the K first by 
calling subroutine RSTRIP1. The subroutine RSTRIP1 solves 
equation 27 for K . RSTRIP1 uses standard technique to find zeros 
of a polynomial. 


For the pulse distribution, the main program solves the matrix 
equation given in equation 36. This is done by computing mutual 
impedances by calling subroutine SRMN. The mutual impedance 
matrix is inverted using MATINV. I are found using MATROW. 

n 

Equation 28 is then solved for K by using subroutine RSTRIP1. 

The characteristic impedance for either distribution is computed 
by calling subroutine MZO. 
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PROGRAM MSTRIP 




C** 

C** 

c** 

c** 

c** 

c** 

c** 

c** 

c** 

c** 

c** 

c** 

c** 

c** 

c mm 

c** 


* * * 

This program calculates characteristic Impedance of a *** 

mlcrostrlp line embedded in a dielectric slab. *** 

**« 


This program gives characteristic impedance of a micro- 
strip line for three types of transverse current distri- 
butions. (1) Uniform, (2) , Nonuniform with singularities 
at the strip edges, (3) General distribution with pulse 
approximations. 


»** 
*** 
* ** 
*** 
*#* 
*** 


This program gives characteristic Impedance of a micro- *** 
strip line as a function of strip width. It also gives *** 

the transverse current distribution. *** 

*** 


DIMENSION BTM (10), CZPAT ( 242 ) , CZMN ( 242 , 242 ) , C I ( 242 ) 

1 , CZT 1 ( 242 , 242 ) , CZT2 (242,242), CZT3 ( 242 , 242 ) , CVT ( 242 ) 

2, CVT1 (242 ) , Cl 1 (242) , CI2(242 ) 

COMPLEX CJ, CJE1 , CZPAT, CZMN, ZPATR, CF1 , CZT1 . CZT2, CZT3, CVT, CVT1 
COMMON CJ, CJE1 , PI , TWOPI , PI02, DELZ, DELD, DELDZ, DELW, XX, B, E, El 


C** 

C** 

C** 

C** 

C** 

C** 

C»* 

C** 

c** 

c** 

c** 

c** 

c»» 

c*» 

c** 

c** 

c** 

c** 

c*» 

c** 

c** 

c** 

c** 

c** 


INPUT DATA 


*••**•••*«**••*••****•••***•** 


NU = INTEGER AND SELECTS TVPE OF CURRENT DISTRIBUTION IN 
TRANSVERSE DIRECTION. 

NU = 1 FOR PULSE DISTRIBUTION 

= 2 FOR UNIFROM DISTRIBUTION 

* 3 FOR NON-UNIFORM DISTRIBUTION 

> or = 4 EXIT 

MM = NUMBER OF PULSES FOR NU = 1 

= ANY NUMBER WHEN NU IS NOT EQUAL TO ONE. 

MAXIMUM VALUE OF MM MUST BE LESS THAN 242. 

ED = DIELECTRIC SLAB THICKNESS IN WAVELENGTH 

EZP = HIGHT OF MICROSTRIP FROM GROUND IN WAVELENGTH 
(EZP MUST BE LESS THAN ED) 


C** E = DIELECTRIC CONSTANT OF SLAB. 

C* * * 

C** WSTRIP = WIDTH OF MICROSTRIP IN WAVELENGTH * 

C** * 

C** WMAX * MAXIMUM WIDTH OF STRIP IN WAVELENGTH * 

C** * 

0****** ********** a******************** *•*******• ******* **»»*****•*• 
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WRITE(5,*) ’GIVE THE VALUE OF NU’ 

WRITE(5, * ) ’ FOR PULSE DISTRIBUTION NU=1, FOR UNIFORM DISTRIBUTION 
1NU=2, FOR NONUNIFORM DISTRIBUTION NU=3, NU GREATER THAN 3 EXIT. ’ 
READ (5, * )NU 

WRITE(5, *)’GIVE THE VALUE OF MM’ 

WRITE(5, * ) ’ MM = 1 TO 242’ 

READ (5, • )MM 

WRITE(5, * ) ’ GIVE THE VALUE OF DILECTRIC SLAB THICKNESS 
1 IN WAVELENGTH' -. A 

READ (5, * )ED 

WRITE(5, * ) ’ GIVE THE POSTION OF STRIP FROM GROUND IN 
1 WAVELENGTH’ 

READ (5, * )EZP 

WRITE(5, * ) ’ GIVE DIELECTRIC CONSTANT AND LOSS TANGENT 
1 OF SLAB MATERIAL’ 

READ(5, * )E, El 

WRITE (5, * ) ’ GIVE STRIP WIDTH, MAXIMUM STRIP WIDTH IN 
1 WAVELENGTH’ 

READ(5, * JWSTRIP, WMAX 
IF(NU. GT. 3)GO TO 105 
DELTX=0 . 001 
CJ=(0. , 1. ) 

PI=2. *ASIN( 1.0) 

PI02=PI*0. 5 
TWOP I =2 . *PI 

CJE1=CJ* (E-CJ*E*E1-1 . ) 

SQRE=SQRT ( E ) 

WRITE! 13, •)’ CHARACTERISTIC IMPEDANCE, PROPAGATION CONSTANT’ 

WRITE! 13, * ) ’ OF MICROSTRIP LINE 

WRITE! 13, TSLAB THICKNESS ED ■ ’ , ED 

WRITE! 13, *)’ POSITION OF STRIP EZP =’,EZP 

IF(NU. EQ. 1 ) WRITE! 13, *)’ DISTRIBUTION PULSE WITH MM =’,MM 

IF(NU. EQ. 2)WRITE( 13, * ) ’ UNIFORM DISTRIBUTION’ 

I F ( NU . EQ . 3 ) WR I TE ( 1 3 , * ) ’ NONUNIFORM D I STR I BUT I ON ’ 

WRITE! 13, * ) ’ **STRIPWIDTH** ’ , ’ ** PROPAGATION CONST. **’ , ’ **IMP. **’ 
103 DELD=TWOPI *ED 
DELZ=TWOPI *EZP 
DELW=TWOP I * WSTR I P 
DELDZ=DELD-DELZ 
DELWY =DELW/MM 
B=0 


*• **** 

C** If dielectric constant of the slab = 1 **** 

C** then avoid surface wave calculation. **** 

q;** •*** 

C* ************************************************************* 

IF(E. EQ. 1. )G0 TO 501 

NE=1 

NNMAX=1 

CALL ZTEMP(BTM, NE, NNMAX, E, DELD) 

B=BTM( 1 ) 
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£***********************************•********♦***************** 

Q*m *** 

C** If NU * 1 calculate mutual impedance matrix using *** 

C** subroutine SRMN. When NU not equal to 1 mutual *** 

C** matrix Is not required. *** 

C** *** 

IFCNU.NE. 1 )G0 TO 502 
501 CALL SRMN (MM, CZPAT ) 

DO 10 1=1, MM 
DO 10 J=1 , MM 
IF( J. LT. I )G0 TO 10 
CZMN(I, J)=CZPAT(J-I+1) 

10 CONTINUE 

DO 20 1=1, MM 
DO 20 J=1 , MM 
IF( J. LE. I)G0 TO 20 
CZMN(J, I ) =CZMN ( I , J ) 

20 CONTINUE 

DO 30 1=1, MM 
DO 30 J = 1 , MM 
CZT1 ( I , J ) =CZMN ( I , J ) 

30 CONTINUE 

EDK=SQRE*DELZ 

DKK=SQRE*DELD 

CF 1 =2 . * C J *S I N ( EDK ) * CEXP ( C J * DELZ ) / ( SQRE * COS ( DKK ) +C J *S I N ( DKK ) ) 

DO 50 1=1, MM 
CVT( I )=-CFl 
50 CONTINUE 

CALL MATINV(CZT1 , MM, CZT2) 

CALL MATRON ( CZT2 , MM, CVT, CVT1 ) 

S1=0. 

DO 301 1=1, MM 
FACT=CABS(CVT1 ( I ) ) 

S1=S1+FACT 

301 CONTINUE 

DO 302 1=1, MM 
CVT1 ( I )=CVT1 ( I )/Sl 

302 CONTINUE 

DO 303 1=1, MM 

Cl ( I )=CABS(CVT1 ( I ) ) 

303 CONTINUE 

£* ******************************************************************** 
£* • * • * 

C** Cl (I) is the current distribution In the transverse direction*** 

c *. *** 

£•************************•***********#******••*********•************* 

WRITE(10,*) ’TRANSVERSE CURRENT DISTRIBUTION’ 

WRITE(10,*) ’****•*• Y/W *****•’, ’*******CURRENT AMP.***** 

DO 304 1=1, MM 
DWY= (2*1-1 ) 

DWY =DW Y * DELW Y * 0 . 5/DELW 
WRITE( 10, * )DWY, Cl (I ) 

304 CONTINUE 
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* 


C***' 

c** 

c** 

c** 

c** 

c** 

c** 

C*** 

502 

C"** 

C** 

C** 

C** 

c** 

c** 

c** 

c** 

c*** 


105 

104 


• •■•••••••••••••******«************* , *********************** l,,ft 

*•* 

The transverse current distribution calculation for. NU =1 *** 

Is completed here. The Cl ( I ) as a function of DWY Is 
outputted at F0R010.DAT file. The Cl (I) Is not computed *** 

when NU Is not equal to one. 


NNMAX=1 

**************************** ********************** ** 

♦ *** 

For the current distribution Cl (I) the SUBROUTINE RSTRIP1 **** 

calculates dominant mode propagation constant for the micro- 

**«* 

strip line. 

RSTRIP1 solves equation (17) for Ke. 

RSTRIP1 must be linked to RTSRIP SUBROUTINUE •*** 

#**» 

«»■••••••••••*••*«•*•••**•************************************* 

CALL RSTRIP1 (NU, NNMAX.MM, DELWY, DELD, DELDZ, DELZ, DELW, 

IE, El, BTM.CI) 

XX=BTM ( 1 ) 

CALL MZO ( NU , DELTX , ZPATR , C I , MM ) 

WRITE( 13, * )WSTRIP, XX, ZPATR 
WSTR I P=WSTR I P+0 . 005 
IF(WSTRIP. GT. WMAX)G0 TO 104 
GO TO 103 

WR I TE ( 5 , * ) ’ CURRENT D I STR I BUT I ON NOT SELECTED PROPERLY ’ 

STOP 

END 
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SUBROUT I NE ZTEMP ( BTM , NE , M AX , E , DK ) 
DIMENSION BTM(l) 

SQRE=SQRT(E) 

DELTA=0. 1 
Bl=l . 0 

B2=l. O+DELTA 

IF(B2. GT. SQRE ) B2=SQRE 

BT=B2 

N=0 .-.js 

10 D1=SQRT(E-B1*B1) 

DK1=DK*D1 
S1=SIN(DK1 ) 

C1=C0S(DK1 ) 

F1=E*SQRT(B1*B1-1 . 0)*C1-D1*S1 
D2=SQRT(ABS(E-B2*B2) ) 

DK2=DK*D2 
S2=SIN (DK2 ) 

C2=C0S ( DK2 ) 

F2=E*SQRT(B2*B2-1. 0)*C2-D2*S2 
FT=F2 

IF(F1. LT. 0. 0. AND. F2. LT. 0. 0) GO TO 15 
I F ( F 1 . GT. 0. 0. AND. F2. GT. 0. 0) GO TO 15 

12 AF1=ABS(F1 ) 

AF2=ABS(F2) 

BP=B1+ (B2-B1 )*AF1/(AF1+AF2) 
IF(ABS(B2-B1 ) . LE. 1 . E-04 ) GO TO 14 
DP=SQRT ( E - BP * BP ) 

DKP=DK*DP 

SP=SIN(DKP) 

CP=COS ( DKP ) 

FP=E*SQRT ( BP*BP-1 . 0 ) *CP-DP*SP 

IF(ABS(FP). LE. 1. E-04) GO TO 14 

IFIFI.LT. 0.0. AND.FP.LT.O.O) GO TO 13 

IF(F1.GT.0.0. AND.FP.GT.O.O) GO TO 13 

B2=BP 

F2=FP 

GO TO 12 

13 B1=BP 
F1=FP 

GO TO 12 

14 N=N+1 
IF(N.GT.MAX) GO TO 17 
BTM(N)=BP 

15 B1=BT 

IF(ABS(B1-SQRE). LE. 1. E-05) GO TO 16 

B2=B1+DELTA 

I F ( B2 . GT . SQRE ) 82=SQRE 

BT=B2 

GO TO 10 

16 NE=N 

17 RETURN 
END 
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SUBROUT I NE SRHN ( MM , ZPATR ) 

DIMENSION U1(3),U2(10),R1(3),R2(10),U(13),R(13) 

D I MENS I ON BTM ( 1 ) , SUMRR ( 242 ) , ZPATR ( 242 ) 

EQUIVALENCE (U1 ( 1 ) , U( 1 ) ) , (U2( 1 ), U(4) ) , (R1 ( 1 ) , R( 1 ) ) . (R2( 1 ) , R(4) ) 
COMMON CJ, CJE1 , PI , TWOPI , PI02, DELZ, DELD, DELDZ, DELW, XX, B, E, El 
COMPLEX CJ, CJE1 , CF1 , SUMRR, ZPATR, CCFA 

DATA Ul/O. 11270166537925, . 5, 0. 88729833462074/, U2/. 01304673574141 

A, . 06746831665550, . 16029521585048, . 28330230293537, . 42556283050918 

B, . 57443716949081, . 71669769706462, . 83970478414951, . 93253168334449 

C, . 98695326425858/, Rl/. 27777777777777, . 44444444444444, . 2777777777 
D7777/, R2/. 03333567215434, . 07472567457529, . 10954318125799, . 134633 
E35965499, . 14776211235737, . 14776211235737, . 13463335965499, . 109543 
F18125799, . 07472567457529, . 03333567215434/ 

DELWY=DELW/MM 
SQRE=SQRT ( E ) 

DWK02=0. 5* DELW Y 
C0NST=377. 7 * DELZ * DELWY/P I 
DO 201 1=1, MM 
SUMRR ( I ) = ( 0 . ,0. ) 

ZPATR ( I ) = (0. ,0. ) 

201 CONTINUE 
NF=1 
NQ=5 

XX1=0.0001 
XX2=1 . -0. 0001 

100 DELT= (XX2 -XXI) /FLOAT (NQ) 

DO 80 K= 1 , NQ 
XI=K-1 

FF =XX 1 +X I * DELT 
DO 80 JL=4, 13 
UU=U( JL)*DELT+FF 
BSQ=UU*UU 
BET=UU 
BETA=UU 

DWKB2=BET A • DWK02 
AKX 1 =BSQ 

IF(E. EQ. 1. )G0 TO 501 
IF(AKX1-1 . ) 10, 10, 12 
10 B1=SQRT ( 1 . -AKX1 ) 

BE=SQRT(E-AKX1 ) 

AZ=BE*DELZ 
AD=DELD*BE 
ADZ=DELDZ*BE 
SINZ=SIN (AZ) 

COSD=COS ( AD ) 

S I ND=S I N ( AD ) 

COSDZ=COS ( ADZ ) 

SINDZ=SIN( ADZ) 

SINZ1=1. 0 

I F ( ABS ( AZ ) . LT . 1 . E-05 ) GO TO 21 
SINZ1=SINZ/AZ 

21 CF1=SINZ1*(BE*C0SDZ+CJ*B1*SINDZ)/(BE*C0SD+CJ*B1*SIND) 

GO TO 20 

12 IF(AKX1-E)14, 14, 15 
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14 B1=SQRT(AKX1-1 . ) 

BE=SQRT ( E-AKX 1 ) 

AZ=BE*DELZ 

AD=BE*DELD 

ADZ=BE*DELDZ 

SINZ=SIN(AZ) 

COSD=COS ( AD ) 

SIND=SIN(AD) 

COSDZ=COS ( ADZ ) 

SINDZ=SIN(ADZ) 

SINZ1=1. 

IF(ABS(AZ) . LT. 1 . E-05)G0 TO 22 
SINZ1=SINZ/AZ 

22 CF1=SINZ1* (BE*COSDZ+Bl *SINDZ)/(BE*COSD+Bl *SIND) 
GO TO 20 

15 B1=SQRT(AKX1-1. ) 

BE=SQRT ( AKX 1 -E ) 

AZ=DELZ*BE 
AD=DELD*BE 
ADZ=DELDZ*BE 
EP1=0. 0 

IF(ABS(AZ).GT20. )G0 TO 23 
EPl=EXP(-2. *AZ) 

23 EPD=0. 

I F ( ABS ( AD ) . GT . 20 . )G0 TO 24 
EPD=EXP(-2. *AD) 

24 EPDZ=0 . 

IF(ABS(ADZ) . GT. 20. }G0 TO 25 
EPDZ=EXP(-2. *ADZ) 

25 CF1=( 1 . -EP1 ) * (BE* ( 1 . +EPDZ) +B1 * ( 1 . -EPDZ) ) 
CF1=CF1/(BE* ( 1 . +EPD)+B1*(1. -EPD) ) 

CFl=CFl/(2. *AZ) 

GO TO 20 

501 IF(AKX1-1. 1502,502,503 

502 BE=SQRT( 1 . -AKX1 ) 

AZ=BE*DELZ 

SINZ=SIN(AZ) 

1F(ABS(AZ) . LT. 1 . E-05)G0 TO 504 
SINZ1=SINZ/AZ 

504 CF 1 =S INZ 1 * ( COS ( AZ ) -C J * S INZ ) 

GO TO 20 

503 BE=SQRT(AKX1-1. ) 

AZ=BE*DELZ 

EZ=0. 

I F ( ABS ( AZ ) . GT . 20 . )G0 TO 505 
EZ=EXP(-2. *AZ ) 

505 CF1=( 1 . -EZ)/(2. *AZ) 

20 DWKA2=DWKB2 

Fl=l . 0 

1 F ( ABS ( DWKA2 ) . LT. 1 . E-05)GO TO 301 
F 1 =S I N ( DWKA2 ) /DWKA2 
301 CCFA=-CJ*C0NST*CF1 
CCFA=CCFA* F 1 * F 1 
DO 202 1=1, MM 
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SUMRR ( I ) =SUMRR ( I ) +CCFA*R( Jl. ) *COS (2. *DWKA2« ( 1 - 1 ) ) 

202 CONTINUE 
80 CONTINUE 

DO 203 1=1, MM 

SUMRR ( I )=SUMRR( I ) *DELT 

203 CONTINUE 

DO 204 1=1, MM 

ZPATR ( I ) =ZPATR ( I ) +SUMRR ( I ) 

204 CONTINUE 

DO 205 1=1, MM 
SUMRR ( I ) = ( 0 . ,0. ) 

205 CONTINUE 
CPOLE=0 . 0001 
NF=NF+1 

IF(NF. EQ. 2)XX1=1. +CPOLE 
IF(NF. EQ . 2 ) XX2=SQRE-CP0LE 
IF(NF. EQ. 2)G0 TO 100 
IF(NF. EQ. 3 ) XX 1 =SQRE+CPOLE 
IF(NF. EQ. 3)XX2=XX1 +1 . 

IF(NF. EQ. 3)G0 TO 100 

IF (XX2. LT. 5. 0. AND. XX2. GT. 0. )DINC=1 . 0 

IF(XX2. LT. 25. 0. AND. XX2.GT. 5. 0)DINC=5. 

IF(XX2. LT. 100. 0. AND. XX2. GT. 25. 0)DINC=25. 0 

IF(XX2. GT. 100. 0)DINC=100. 0 

XX1=XX2 

XX2=XX1+DINC 

IF(XX2. LT. 2000. 0)GO TO 100 

RETURN 

END 
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SUBROUT I NE MAT I NV ( SUM 1 , N , SUM2 ) 

DIMENSION I ST (242) , JST(242 ) 

COMPLEX SUM1 ( 242, 242 ) , SUM2 ( 242, 242 ) 

DO 1 1=1, N 
IST( I )=0 
1 JST ( I ) =0 

DO 30 IND=1 , N 
H=0. 0 

DO 10 1=1, N 
DO 3 KCH=1 , N 

IF t IST(KCH) . EQ. I )G0 TO 10 

3 CONTINUE 
DO 9 J=1,N 
DO 4 KCH=1 , N 

IF( JST(KCH) . EQ. J ) GO TO 9 

4 CONTINUE 

IF(H. GE. CABS(SUM1 ( I , J) ) ) GO TO 9 
H=CABS(SUM1 ( I , J ) ) 

IM=I 

JM=J 

9 CONTINUE 
10 CONTINUE 

IST(IND)=IM 
JST( IND)=JM 
DO 20 1=1, N 
DO 20 J=1 , N 
IF ( I-IM)5, 6, 5 

6 SUM2(I, J)=-SUM1(I, J)/SUM1(IM, JM) 

GO TO 20 

5 IF( J-JM)7, 8, 7 

8 SUM2 ( I , J)=SUM1 ( I , J )/SUMl ( IM, JM ) 

GO TO 20 

7 SUM2( I , J )=SUM1 ( I , J )-SUMl ( I , JM)*SUM1 ( IM, J )/SUMl ( IM, JM) 
20 CONTINUE 

SUM2( IM, JM)=1 . /SUM1 ( IM, JM) 

DO 25 1=1, N 
DO 25 J=1 , N 

25 SUM1 ( I , J ) =SUM2 ( I , J ) 

30 CONTINUE 

DO 35 1=1, N 

I F ( IST( I ) . EQ. JST( I ) )G0 TO 35 
DO 45 J-l.N 

45 SUM2( JST( I ) , J)=SUM1 ( IST( I ) , J ) 

35 CONTINUE 

DO 26 1=1, N 
DO 26 J=1 , N 

26 SUM 1 ( I , J ) =SUM2 ( I , J ) 

DO 50 1=1, N 

I F ( IST( I ) . EQ. JST ( I ) )G0 TO 50 
DO 55 J=1 , N 

55 SUM2(J, IST(I))=SUM1(J, JST(I)) 

50 CONTINUE 
RETURN 
END 
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SUBROUT I NE MATROW ( SUM 2 , N , SUM, SUM5 ) 
COMPLEX SUM2 ( 242 , 242 ) , SUM ( 242 ) , SUMS ( 242 ) 
DO 110 1=1, N 
SUMS ( I ) = ( 0 . ,0. ) 

110 CONTINUE 

DO 120 1=1, N 
DO 120 J=1 , N 

SUM5 ( I )=SUM2 ( I , J ) *SUM( J )+SUM5 ( I ) 

120 CONTINUE v , 

RETURN 

END 
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SUBROUTINE RSTRIP1 (NU, NNMAX, MM, DELWY, DK, DZP, ZPK, DWSTRIP, 

IE, El , BTM, Cl ) 

DIMENSION BTM( 1 ) , Cl (242 ) 

COMPLEX ZPATRR 
SQRE=SQRT(E) 

N=0 

DELTA=0. 1 
B1=SQRE 

B2=B1 -DELTA .,i 

10 CALL RSTRIP(NU, MM, DELWY.DK, DZP, ZPK, DWSTRIP.E, El, Bl, ZPATRR, Cl) 

F1=REAL( ZPATRR) 

CALL RSTRIP(NU, M, DELWY, DK, DZP, ZPK, DWSTRIP, E, El , B2, ZPATRR, Cl ) 
F2=REAL( ZPATRR) 

FT=F2 

IF(F1. LT.O.O. AND. F2.LT.0.0) GO TO 15 
IF(F1.GT. 0. 0. AND.F2.GT. 0.0) GO TO 15 

12 AF1=ABS(F1 ) 

AF2=ABS(F2) 

BP=B1-(B1-B2)*AF1/(AF1+AF2) 

IF(ABS(B2-B1).LE. l.E-04) GO TO 14 

CALL RSTRIPCNU.MM, DELWY, DK, DZP, ZPK, DWSTRIP, E, El , BP, ZPATRR, Cl ) 
FP=REAL( ZPATRR) 

I F ( ABS ( FP ) . LE . l.E-04) GO TO 14 

IF (FI . LT.O.O. AND.FP.LT. 0.0) GO TO 13 

IF(Fl.GT.O. 0. AND.FP.GT. 0.0) GO TO 13 

B2=BP 

F2=FP 

GO TO 12 

13 B1=BP 
F1=FP 

GO TO 12 

14 IF(ABS(FP) . LT. 1 . E-04)G0 TO 19 
GO TO 15 

19 N=N+1 

BTM(N)=BP 
IFCN.GE. 1 )G0 TO 17 

15 B1=B2 
B2=B1 -DELTA 
GO TO 10 

16 NE=N 

17 RETURN 
END 
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SUBROUTINE RSTRIP(NU, MM, DELWY, DK, DZP, ZPK, DWSTRIP, E, El , 

1 XX, ZPATR, Cl ) 

DIMENSION U1 (3),U2(10),R1(3),R2(10),U(13),R(13), BTM( 1 ) , Cl (242 ) 
EQUIVALENCE (U1 ( 1 ) , U( 1 ) ) , (U2 ( 1 ) , U(4) ) , (R1 ( 1 ) , R( 1 ) ) , (R2( 1 ) , R(4) ) 
COMPLEX CJ, CJE1 , CF1 , CF2, SUMRR, ZPATR, CCFA 

DATA Ul/O. 11270166537925, . 5, 0. 88729833462074/, U2/. 01304673574141 

A, . 06746831665550, . 16029521585048, . 28330230293537, . 42556283050918 

B, . 57443716949081 , . 71669769706462, . 83970478414951 , . 93253168334449 

C, . 98695326425858/, Rl/. 27777777777777, . 44444444444444, . 2777777777 
D7777/, R2/. 03333567215434, . 07472567457529, . 10954318125799, . 134633 
E35965499, . 14776211235737, . 14776211235737, . 13463335965499, . 109543 
F18125799, . 07472567457529, . 03333567215434/ 

CJ=(0. , 1. ) 

PI =2. ■ASINU. 0) 

TWOP I =2 , "PI 
PI02=PI/2. 

SQRE=SQRT ( E ) 

DWK02=0. 5 "DWSTRIP 
DWKY2=0. 5 "DELWY 
CJE1=CJ* (E-CJ*E1 "E-l . ) 

SUMRR= ( 0 . ,0. ) 

ZPATR=(0. ,0. ) 

CONTINUE 

NF=1 

NQ=5 

XX1=0.0001 
XX2=1 . -0. 0001 

100 DELT=(XX2-XX1 ) /FLOAT (NQ) 

DO 80 K=1 , NQ 
XI=K-1 

FF=XX1+XI "DELT 
DO 80 JL=4, 13 
UU=U ( JL ) "DELT +FF 
BSQ=UU*UU 
BET=UU 
BETA=UU 

DWKB2=BETA*DWK02 
DYKB2=BETA*DWKY2 
AKX1=XX*XX+BSQ 
IF(AKX1-1 . ) 10, 10, 12 
10 B1=SQRT( 1 . -AKX1 ) 

BE=SQRT(E-AKX1 ) 

AZ=BE*ZPK 

AD=DK*BE 

ADZ=DZP"BE 

SINZ=SIN(AZ) 

COSD=COS ( AD ) 

S I ND=S I N ( AD ) 

COSDZ=COS ( ADZ ) 

SINDZ=SIN ( ADZ) 

SINZ1=1 . 0 

IF(ABS(AZ). LT. 1. E-05)G0 TO 21 
SINZ1=SINZ/AZ 

21 CF1=SINZ1 " (BE*C0SDZ+CJ*B1 *SINDZ)/(BE*C0SD+CJ*B1 "SIND) 
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CF2=SINZ1*SINZ*BE*BE/(BE*C0SD+CJ*B1*SIND) 
CF2=CF2/ ( E*B1 ‘COSD+CJ *BE*S I ND ) 

GO TO 20 

12 IF(AKX1-E}14, 14. 15 

14 B1=SQRT(AKX1 -1 . ) 

BE=SQRT(E-AKX1 ) 

AZ=BE*ZPK 

AD=BE*DK 

ADZ=BE*DZP 

SINZ=SIN(AZ) 

COSD=COS ( AD ) 

S I ND=S I N ( AD ) 

COSDZ=COS ( ADZ ) 

SINDZ=SIN(ADZ) 

SINZ1=1 . 

IF(ABS(AZ). LT. 1. E-05)G0 TO 22 
SINZ1=SINZ/AZ 

22 CF1=SINZ1* (BE*C0SDZ+B1 *SINDZ)/(BE*C0SD+B1 *SIND) 
CF2=SINZ1*SINZ*BE*BE/(BE*C0SD+B1 *COSD) 
CF2=CJ*CF2/ (B1*E*C0SD-BE*SIND ) 

GO TO 20 

15 B1=SQRT(AKX1-1. ) 

BE=SQRT ( AKX 1 -E ) 

AZ=ZPK*BE 
AD=DK*BE 
ADZ=DZP*BE 
EP1=0. 0 

IF(ABS(AZ). GT. 20. )G0 TO 23 
EPl=EXP(-2. *AZ) 

23 EPD=0 . 

IF(ABS(AD) . GT. 20. )G0 TO 24 
EPD=EXP(-2. *AD) 

24 EPDZ=0 . 

IF(ABS(ADZ) . GT. 20. )G0 TO 25 
EPDZ=EXP(-2. *ADZ) 

25 CF1 = ( 1 . -EP1 ) * (BE* ( 1 . +EPDZ)+B1 * ( 1 . -EPDZ) ) 
CF1=CF1/(BE* ( 1 . +EPD) +B1* ( 1 . -EPD) ) 

CFl=CFl/(2. *AZ) 

CF2=-CJ*BE*BE* ( 1 . -EP1 ) * ( 1 . -EP1 ) 'EPDZ/AZ 
CF2=CF2/(BE* [ 1 . +EPD )+Bl* ( 1 . -EPD) ) 

CF2=CF2/(E*B1* ( 1 . +EPD)+BE* ( 1 . -EPD) ) 

20 DWKA2=DWKB2 

IF(NU. NE. 1 )G0 TO 501 
FI 1=SIN (DYKB2 ) 

FI = 1 . 

IF(ABS(DWKA2). LT. 1. E-05)G0 TO 301 
F 1 =S I N ( DWKA2 ) /DWKA2 

301 F2=l . 

IF(ABS(DYKB2) . LT. 1. E-05)GO TO 302 
F2=F1 1/DYKB2 

302 GO TO 508 

501 F3=0. 5*DWKA2 

IF(ABS(F3) . LT. 1 . E-05)GO TO 506 
F3=SIN(F3)/F3 
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F4=l . 

I F ( ABS ( DWKA2 ) . LT. 1 . E-05)G0 TO 507 
F 4=S I N ( DWKA2 ) /DWKA2 

507 Fl=l . 

F2=l . 

IF(ABS(DWKA2 ) . LT. l.E-05)G0 TO 503 
IF(NU. EQ. 2)F2=F4 

IF(NU. EQ. 3 ) CALL BJ0R(DWKA2, F2, I ERR) 

IF(NU. EQ. 4)F2=2. *F4-F3*F3 •-/. 

503 F2=F2*F2 

508 CCFA=( (E-XX*XX) *CF1-CJE1 *XX*XX*CF2) 
CCFA=CCFA*F1 *F2 

IFCNU. NE. 1 )G0 TO 504 
S1=0. 

DO 303 1=1, MM 
XARG=MM/2 
XARG1=I 
XARG2=0. 5 

XARG=XARG- ( XARG1 -XARG2 ) 

FACT=CI ( I ) *C0S(2. *DYKB2*XARG) 

S1=S1+FACT 
303 CONTINUE 
GO TO 505 

504 Sl-1. 

505 SUMRR=SUMRR + CCF A • R ( J L ) * S 1 
80 CONTINUE 

SUMRR=SUMRR * DELT 
ZPATR=ZPATR+SUMRR 
SUMRR= ( 0 . ,0. ) 

CP0LE=0 . 0001 
NF=NF+1 

IF(NF. EQ. 2)XX1 = 1 . +CPOLE 
IF(NF. EQ. 2 ) XX2=SQRE-CP0LE 
IF(NF. EQ. 2)G0 TO 100 
IF(NF. EQ. 3 ) XXI =SQRE+CPOLE 
IF(NF. EQ. 3)XX2=XX1+1. 

IF(NF. EQ. 3)G0 TO 100 

IF(XX2. LT. 5. 0. AND. XX2. GT. 0. )DINC=1. 0 

IF(XX2. LT. 25. 0. AND. XX2. GT. 5. 0)DINC=5. 

IF(XX2. LT. 100. 0. AND. XX2. GT. 25. 0)DINC=25. 0 

IF(XX2. GT. 100. 0)DINC=100. 0 

XX1=XX2 

XX2=XX1+DINC 

IF(XX2. LT. 1000. 0)GO TO 100 

RETURN 

END 
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SUBROUT I NE MZO ( NU , DELTX , ZP ATR , C I , MM ) 

C** «.**•*•*•*•******•* ****************** 

C** This program computes the characteristic Impedance of micro 

C** line. This program calls subroutine DFKX to calculate 
C** derivatives of the integrand appearing in equations (12), 
C** (13), and (14). 

C** 

C** 

C** INPUT PARAMETERS: v , 

C** 

C** NU = 1 FOR PULSE DISTRIBUTION 
C** = 2 FOR UNIFORM DISTRIBUTION 

C** =3 FOR NONUNIFORM DISTRIBUTION 

C** 

C** DELTX = A CONSTANT IN THE RANGE 0.01-0.001 

C** 

c** Cl = CURRENT DISTRIBUTION AMP. 

C** IT IS ZERO IF NU=2 OR NU=3 

C** 

C** MM = NUMBER OF PULES 

C** 

C** OUTPUT PARAMETERS: 


C** 

C** ZPATR = CHARACTERISTIC IMPEDANCE 

C»* 

£***#************** ********* ********** 


OF MICROSTRIP LINE * 

* 


DIMENSION U1 (3),U2(10),R1 ( 3) , R2( 10) , U( 13) , R(13) 

DIMENSION BTM( 1 ) , Cl (242) 

EQUIVALENCE (U1 ( 1 ) , U( 1 ) ) , (U2( 1 ) , U(4 ) ) . (R1 ( 1 ) , R( 1 ) ) , (R2( 1 ) , R(4 ) ) 
COMMON CJ, CJE1 , PI , TWOPI , PI02, DELZ, DELD, DELDZ, DELW, XX, B, E, El 


COMPLEX CJ, CJE1 , CF1 , CF2, SUMRR, ZPATR, CCFA, DF 

DATA Ul/O. 11270166537925, .5, 0.88729833462074/, U2/. 01304673574141 
A, . 06746831665550, . 16029521585048, . 28330230293537, . 42556283050918 
B ; . 57443716949081, . 71669769706462, . 83970478414951, . 93253168334449 
c’ 98695326425858/, Rl/. 27777777777777, . 44444444444444, . 2777777777 
D7777/, R2/. 03333567215434, . 07472567457529, . 10954318125799, . 134633 
E35965499, . 14776211235737, . 14776211235737, . 13463335965499, , 109543 
F18125799, . 07472567457529, . 03333567215434/ 

CPOLE=0 . 0001 
SQRE=SQRT ( E ) 

DELWY=DELW/MM 
DWK02=0. 5* DELW 
DYK02=0. 5*DELWY 
CONST=-377. 7*DELZ/(PI*E) 


SUMRR=(0. ,0. ) 
ZPATR= ( 0 . ,0. ) 


NF=1 


NQ=5 

XXI =0.0001 
XX2=1 . -0.0001 

100 DELT=(XX2-XX1 ) /FLOAT (NQ) 
DO 80 K=1 , NQ 
XI=K-1 

FF=XX1+XI*DELT 
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DO 80 JL=4, 13 
UU=U ( JL ) * DELT +FF 
BSQ=UU*UU 
BET=UU 
BETA=UU 

DWKB2=BET A * DWK02 
DYKB2=BETA*DYK02 

CALL DFKX(XX, DELTX, DELD. DELDZ, DELZ, E. CJ. CJE1 , BETA, D. ) 
AKX1=XX*XX+BSQ •.,> 

IF(AKX1-1. )10, 10, 12 
10 B1=SQRT( 1 . -AKX1 ) 

BE=SQRT ( E-AKX1 ) 

AD=DELD*BE 
AZ=DELZ*BE 
ADZ=DELDZ*BE 
SINZ=SIN ( AZ) 

SIND=SIN(AD) 

COSD=COS ( AD ) 

SINDZ=SIN (ADZ) 

COSDZ=COS ( ADZ ) 

SINZ1=1 . 

IF( ABS( AZ) . LT. 1 . E-05 )G0 TO 201 
S I NZ 1 =S I NZ/ AZ 

201 CF1=SINZ1 * (BE*C0SDZ+CJ*B1 *SINDZ) 
CF1=CF1/(BE*C0SD+CJ*B1 *S JND) 

CF2=BE*BE*SINZ1 *SINZ/ (BE*C0SD+CJ*B1 *SIND) 

CF2=CF2/ ( E* B 1 * COSD+C J * BE • S I ND ) 

GO TO 20 

12 IF(AKX1-E)14, 14, 18 

14 B1=SQRT(AKX1-1 . ) 

BE=SQRT ( E -AKX1 ) 

AD=DELD*BE 
AZ=DELZ*BE 
ADZ=DELDZ*BE 
SINZ=SIN ( AZ) 

SIND=SIN( AD) 

COSD=COS ( AD ) 

S I NDZ=S I N ( ADZ ) 

COSDZ=COS ( ADZ ) 

SINZ1=1 . 

IF( ABS( AZ) . LT. 1 . E-05 ) GO TO 202 
SINZ1=SINZ/AZ 

202 CF 1 =S I NZ 1 * ( BE * COSDZ+B 1 * S I NDZ ) / ( BE* COSD+B 1 * S I ND ) 
CF2=CJ*BE*BE*SINZ1 *SINZ/ (BE*C0SD+B1 *SIND) 

CF2=CF2/ (E*B1 *COSD-BE*SIND) 

GO TO 20 

18 B 1 =SQRT ( AKX 1 - 1 . ) 

BE=SQRT ( AKX 1 -E ) 

AD=DELD*BE 
AZ=DELZ*BE 
ADZ=DELDZ*BE 
AZ2=2. *AZ 
AD2=2 . *AD 
ADZ2=2 . *ADZ 
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203 


EZ=0. 

I F ( ABS ( AZ2 ) . GT . 15. )G0 TO 203 
EZ=EXP ( -AZ2 ) 

EDD=0 . 

IF(ABS(AD2).GT. 15. )G0 TO 204 
EDD=EXP( -AD2 ) 

204 EDZ=0 . 

IF(ABS(ADZ2) . GT. 15. )G0 TO 205 
EDZ=EXP ( -ADZ2 ) 

205 CF1 = ( 1 . -EZ)* (BE* ( 1 . +EDZ)+B1*(1. -EDZ) ) 
CF1=CF1/(BE*(1. +EDD)+B1*(1. -EDD) ) 

CF1=CF1/AZ2 

CF2=BE*BE* ( 1 . -EZ) * ( 1 . -EZ ) *EDZ/AZ 
CF2=-CJ*CF2/ (BE* ( 1 • +EDD)+B1*(1. -EDD)) 
CF2=CF2/(E*B1* ( 1 . +EDD ) +BE* ( 1 . -EDD ) ) 

20 IF(NU.NE. 1)G0 TO 501 
S1=0. 

DO 302 1 = 1 , MM 
DO 302 J=1 , MM 
XARG=( I-J ) 

FACT=CI (I )*CI ( J)*C0S(2. *DYKB2*XARG) 
S1=S1+FACT 

302 CONTINUE 
F2-1. 

IF(ABS(DYKB2) . LT. 1 . E-05)G0 TO 303 
F2=S I N ( DYKB2 ) /DYKB2 

303 GO TO 502 

501 Sl=l. 

F3=0. 5*DWKB2 

I F ( ABS ( F3 ) . LT . 1 . E-05 ) GO TO 507 
F3=SIN (F3)/F3 

507 F2=l . 

I F ( ABS ( DWKB2 ) . LT . 1 . E-05 ) GO TO 508 
F2=S I N ( DWKB2 ) /DWKB2 

508 IF(NU. EQ. 2)F2=F2 

IF(NU. EQ. 3)CALL BJ0R(DWKB2, F2, I ERR) 

IF(NU. EQ. 4)F2=2. *F2-F3*F3 

502 CCFA=( (E-XX*XX)*CF1-CJE1*XX*XX*CF2-XX*DF)/XX 
CCFA=CONST *DF*F2*F2*S1 
SUMRR=SUMRR+CCFA*R ( JL ) 

80 CONTINUE 

SUMRR=SUMRR * DELT 
ZPATR=ZPATR+SUMRR 
SUMRR= ( 0 . ,0. ) 

NF-NF+1 

IF(NF. EQ. 2)XX1=1 . +CPOLE 
IF(NF. EQ. 2)XX2=B-CPOLE 
IF(NF. EQ. 2)G0 TO 100 
IF(NF. EQ. 3)XX1=B+CP0LE 
IF(NF. EQ. 3 ) XX2=SQRE-CP0LE 
IF(NF. EQ. 3)G0 TO 100 
IF(NF. EQ. 4 )XX1=SQRE+CP0LE 
IF(NF. EQ. 4)XX2=XX2+1 . 0 
IF(NF. EQ. 4)G0 TO 100 
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IF(XX2. LT. 5. 0. AND. XX2. GT. 0. )DINC=1 . 0 
IF(XX2. LT. 25. 0. AND. XX2. GT. 5. 0)DINC=5. 
IF(XX2. LT. 100. 0. AND. XX2. GT. 25. 0)DINC=25. 0 
IF(XX2. GT. 100. 0)DINC=100. 0 
XX1=XX2 
XX2=XX1 +DINC 

IF(XX2. LT. 3000. 0)G0 TO 100 
ZPATR=ZPATR/4. 

RETURN vt 

END 
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SUBROUTINE DFKX(XX, DELTX, DK, D ZP, ZPK, E, CJ, CJE1 , BETA, DF) 

C * * 

C** 

C** This program calculates the derivative of Integrand 

C** appearing In equations (12), (13), and (14). 

C"" 

c **.*************«**** 

DIMENSION DF1 (2) 

COMPLEX DF.DFl 
COMPLEX CJ, CJE1 , CF1 , CF2 
AKX 1 =XX-DELTX 
1 = 1 

BSQ=BET A * BET A 
100 AB1=BSQ+AKX1 "AKX1 
IF(AB1 -1 . ) 10, 10, 12 
10 B1=SQRT( 1 . -AB1 ) 

BE1=SQRT(E-AB1 ) 

AD=BE1*DK 
AZ=BE1 "ZPK 
ADZ=BE1 "ADZ 
SINZ=SIN(AZ) 

SIND=SIN ( AD) 

C0SD=C0S ( AD ) 

SINDZ=SIN(ADZ) 

COSDZ=COS ( ADZ ) 

SINZ1=1 . 0 

IF(ABS(AZ).LT. l.E-05)G0 TO 21 
SINZ1=SINZ/AZ 

21 CF1=SINZ1 * (BE1"C0SDZ+CJ*B1 *SINDZ)/(BE1*C0SD+CJ*B1 *SIND) 
CF2=SINZ1*SINZ*BE1*BE1/(BE1*C0SD+CJ*B1*SIND) 
CF2=CF2/(E*B1*C0SD+CJ*BE1*SIND) 

GO TO 20 

12 IF(AB1-E)14, 14, 15 

14 B 1 =SQRT ( AB 1 - 1 . ) 

BE1=SQRT(E-AB1 ) 

AD=BE1*DK 
AZ=BE1 "ZPK 
ADZ=BE1*DZP 
SINZ=SIN(AZ) 

SIND=SIN(AD) 

COSD=COS(AD) 

SINDZ=SIN( ADZ) 

COSDZ=COS ( ADZ ) 

SINZ1=1 . 

IF(ABS(AZ). LT. 1. E-05)G0 TO 22 
SINZ1=SINZ/AZ 

22 CF1=SINZ1* (BE1*C0SDZ+B1*SINDZ)/(BE1 *C0SD+B1*SIND) 
CF2=SINZ1 *SINZ*BE1 "BE1/ (BE1 "C0SD+B1 "SIND) 
CF2=CJ*CF2/(E*B1*C0SD-BE1*SIND) 

GO TO 20 

15 B1=SQRT(AB1-1. ) 

BE1=SQRT(AB1-E) 

AZ=BE1 "ZPK 
AD=BE1*DK 
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ADZ=BE1 *DZP 
EP1=0. 0 

IF(ABS(AZ) . GT. 20. )GO TO 23 
EPl=EXP(-2. *AZ ) 

23 EPD=0. 0 

IF(ABS(AD) . GT. 20. )G0 TO 24 
EPD=EXP ( -2 . *AD) 

24 EPDZ=0. 0 

I F ( ABS ( ADZ ) . GT . 20 . )G0 TO 25 
EPDZ=EXP(-2. *ADZ) 

25 CF1=( 1 . -EP1 ) " (BE1 * ( 1 . +EPDZ)+B1 * ( 1 . -EPDZ) ) 
CF1=CF1/(BE1* ( 1 . +EPD)+B1* ( 1 . -EPD) ) 
CFl=CFl/(2. *AZ) 

CF2=-CJ*BE1 *BE1* ( 1 . -EP1 ) * ( 1 . -EP1 )*EPDZ/AZ 
CF2=CF2/(BE1 * ( 1 . +EPD) +B1* ( 1 • -EPD) ) 

CF2=CF2/ ( E* B 1 * ( 1 . +EPD ) +BE1 * ( 1 . -EPD ) ) 

20 DF 1 ( I ) = ( E- AKX 1 * AKX 1 ) * CF 1 -C J E 1 * AKX 1 * AKX 1 *CF2 

AKX1=XX+DELTX 
1=1 + 1 

IF(I.LE. 2)G0 TO 100 
DF=(DF1 (2)-DFl ( 1 ) )/DELTX 
RETURN 
END 
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